% this code is for Jumps and Post-FOMC Announcement Returns in Currency Markets
% this gives data during 12hr windows and data for figure 1 and table 1
% (the figure is constructed via MS Excel)

clear
clc

load(using data after filtering intraday data)

%% 1. 12 hr windows and return patterns

% 1) making fomc / non-fomc windows

fomcwindow=nan(nt,2);
for j=1:nt-1
    if usfomc(j,1)==1 && usfomc(j+1,1)==1
        fomcwindow(j-npday*2+1:j+npday*2,1)=1;
        fomcwindow(j-npday*2+1:j+npday*2,2)=[1:1:4*npday]';
    end
end
clear j
% usfomc: dummy for FOMC announcements

nonfomcwindow=nan(nt,2);
for j=npday*4+1:nt-npday*4-1
    if date(j,11)==1
        if date(j,2)==1930 && nansum(fomcwindow(j:j+npday*4-1,1))==0 && ...
                nansum(nonfomcwindow(j-npday*2:j,1))==0 && ...
                nansum(fomcwindow(j-npday*2:j,1))==0
            nonfomcwindow(j-npday*2:j+npday*2-1,1)=1;
            nonfomcwindow(j-npday*2:j+npday*2-1,2)=[1:1:4*npday]';
        end
    end
    if date(j,11)==0
        if date(j,2)==1830 && nansum(fomcwindow(j:j+npday*4-1,1))==0 && ...
                nansum(nonfomcwindow(j-npday*2:j,1))==0 && ...
                nansum(fomcwindow(j-npday*2:j,1))==0
            nonfomcwindow(j-npday*2:j+npday*2-1,1)=1;
            nonfomcwindow(j-npday*2:j+npday*2-1,2)=[1:1:4*npday]';
        end
    end
end
clear j

% 2) weekly window: date(:,5)
% starting at 1999.01.06 (day 6)

% 3) 12 hr window: numbering

tempa=1;
hlfwinnumber=nan(nt,1);

for j=1:nt
    if (nonfomcwindow(j,1)>0)+(fomcwindow(j,1)>0)>0 && ...
            (nonfomcwindow(j,2)==1)+(fomcwindow(j,2)==1)>0
        hlfwinnumber(j:j+npday/2-1,1)=tempa;
        hlfwinnumber(j+npday/2:j+npday-1,1)=tempa+1;
        hlfwinnumber(j+npday:j+npday*(3/2)-1,1)=tempa+2;
        hlfwinnumber(j+npday*(3/2):j+2*npday-1,1)=tempa+3;
        hlfwinnumber(j+2*npday:j+(5/2)*npday-1,1)=tempa+4;
        hlfwinnumber(j+(5/2)*npday:j+3*npday-1,1)=tempa+5;
        hlfwinnumber(j+3*npday:j+(7/2)*npday-1,1)=tempa+6;
        hlfwinnumber(j+(7/2)*npday:j+4*npday-1,1)=tempa+7;
        tempa=tempa+4*2;
    end        
end
clear j

nhlf=max(hlfwinnumber);

winstart=nan(nhlf,1);
for hw=1:nhlf
    for k=2:nt
        if hlfwinnumber(k,1)==hw && (hlfwinnumber(k-1,1)==hw-1)+...
                ((hlfwinnumber(k-1,1)>=0)+(hlfwinnumber(k-1,1)>=0)==0)==1
            winstart(hw,1)=k;
        end
    end
end
clear hw k tempa temprk

% cumulative excess returns

temp=nan(nt,nc+1);
nontemp=nan(nt,nc+1);

for i=1:nc+1
    temp(:,i)=fomcwindow(:,1);
    nontemp(:,i)=nonfomcwindow(:,1);
end
clear i

fomcrxe=rxe.*rxie.*temp;
nonfomcrxe=rxe.*rxie.*nontemp;

clear temp nontemp

sumfomcwindowret=zeros(npday*4,nc+1);
sumnonwindowret=zeros(npday*4,nc+1);

for j=1:npday*4
    temp=nan(nt,nc+1);
    nontemp=nan(nt,nc+1);
    for k=1:nt
        if fomcwindow(k,2)==j
            temp(k,:)=fomcrxe(k,:);
        end
        if nonfomcwindow(k,2)==j
            nontemp(k,:)=nonfomcrxe(k,:);
        end
    end
    sumfomcwindowret(j,:)=nanmean(temp,1);
    sumnonwindowret(j,:)=nanmean(nontemp,1);
end
clear k j temp nontemp

% sum
sumfomcwindowcumret=zeros(npday*4+1,nc+1);
sumnonwindowcumret=zeros(npday*4+1,nc+1);

for j=2:npday*4+1
    sumfomcwindowcumret(j,:)=sumfomcwindowcumret(j-1,:)...
        +sumfomcwindowret(j-1,:);
    sumnonwindowcumret(j,:)=sumnonwindowcumret(j-1,:)...
        +sumnonwindowret(j-1,:);
end
clear j

% figure 1
repbasefigure=zeros(1+npday*4,2*(nc+1));
for i=1:nc+1
    repbasefigure(:,2*i-1)=sumfomcwindowcumret(:,i);
    repbasefigure(:,2*i)=sumnonwindowcumret(:,i);
end
clear i

%% 2. 12 hr window returns and decomposed volatilities

% 1) 12 hr returns

% making lre without jump & only jump & only jup or jdn
tnjlre=nan(nt,nc+1);
tjjlre=nan(nt,nc+1);
tjuplre=nan(nt,nc+1);
tjdnlre=nan(nt,nc+1);

for j=1:nt
    for i=1:nc+1
        if jp(j,i)<1
            tnjlre(j,i)=lre(j,i);
        end
        if jp(j,i)>0
            tjjlre(j,i)=lre(j,i);
        end
        if jup(j,i)>0
            tjuplre(j,i)=lre(j,i);
        end
        if jdn(j,i)>0
            tjdnlre(j,i)=abs(lre(j,i));
        end
    end
end
clear i j

% making hlfday window based return and moments

hwrxe=nan(nhlf,nc+1);
hwlre=nan(nhlf,nc+1);
hwvol=nan(nhlf,nc+1);
hwskew=nan(nhlf,nc+1);
hwkurt=nan(nhlf,nc+1);
hwjrvol=nan(nhlf,nc+1);
hwjvol=nan(nhlf,nc+1);
hwjupv=nan(nhlf,nc+1);
hwjdnv=nan(nhlf,nc+1);
hwjv=nan(nhlf,nc+1);

for hw=1:nhlf
    tempie=lreie(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    temprxie=rxie(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    temprx=rxe(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    templr=lre(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempnjlr=tnjlre(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempjjlr=tjjlre(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempjuplr=tjuplre(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempjdnlr=tjdnlre(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    
    for i=1:nc+1
        ttempie=tempie(:,i);

        ttemprxie=temprxie(:,i);
        ttemprx=temprx(:,i);
        ttemplr=templr(:,i);
        ttempnjlr=tempnjlr(:,i);
        ttempjjlr=tempjjlr(:,i);
        ttempjuplr=tempjuplr(:,i);
        ttempjdnlr=tempjdnlr(:,i);
        
        ttempie(isnan(ttempie))=[];
        ttemprxie(isnan(ttemprxie))=[];
        ttemprx(isnan(ttemprx))=[];
        ttemplr(isnan(ttemplr))=[];
        ttempnjlr(isnan(ttempnjlr))=[];
        ttempjjlr(isnan(ttempjjlr))=[];
        ttempjuplr(isnan(ttempjuplr))=[];
        ttempjdnlr(isnan(ttempjdnlr))=[];
        
        if size(ttempie,1)>npday/4-1 && size(ttemprx,1)>npday/4-1 && ...
                size(ttemplr,1)>npday/4-1 && size(ttemprxie,1)>npday/4-1
            hwrxe(hw,i)=mean(ttemprx)/100*npday/2;    % percentage per 12hr
            hwlre(hw,i)=mean(ttemplr)*npday/2*100;
            hwvol(hw,i)=sqrt(mean(ttemplr.^2)*npday/2);
            hwskew(hw,i)=(sqrt(npday/2)*mean(ttemplr.^3)*npday/2)/...
                hwvol(hw,i)^3;
            hwkurt(hw,i)=(npday/2*mean(ttemplr.^4)*npday/2)/hwvol(hw,i)^4;
            
            hwjrvol(hw,i)=sqrt(((npday/2)/size(ttempie,1))*...
                sum(ttempnjlr.^2));
            hwjvol(hw,i)=sqrt(sum(ttempjjlr.^2)*(npday/2)/size(ttempie,1));
            hwjupv(hw,i)=sqrt(sum(ttempjuplr.^2)*(npday/2)/size(ttempie,1));
            hwjdnv(hw,i)=sqrt(sum(ttempjdnlr.^2)*(npday/2)/size(ttempie,1));
            hwjv(hw,i)=(sum(ttempjuplr.^2)+sum(ttempjdnlr.^2))*...
                (npday/2)/size(ttempie,1);
            if hwjv(hw,i)>hwvol(hw,i)^2
                hwjv(hw,i)=hwvol(hw,i)^2;
            end
        end
    end
    hw
end
clear hw k i
clear tnjlre tjjlre tjuplre tjdnlre tnjlrem tjjlrem tjuplrem tjdnlrem
clear ttempie ttemprx ttemplr ttempnjlr ttempjjlr ttempjuplr ttempjdnlr
clear ttempnjlrm ttempjjlrm ttempjuplrm ttempjdnlrm
clear tempie temprx templr tempnjlr tempjjlr tempjuplr tempjdnlr
clear tempnjlrm tempjjlrm tempjuplrm tempjdnlrm ttemprxie

hwvol=100*hwvol;
hwjrvol=100*hwjrvol;
hwjupv=100*hwjupv;
hwjdnv=100*hwjdnv;
hwjvol=100*sqrt(hwjv);

% 2) 12 hr residual returns

% making resilr/resirx without jump & only jump & only jup or jdn
tnjlre=nan(nt,nc);
tjjlre=nan(nt,nc);
tjuplre=nan(nt,nc);
tjdnlre=nan(nt,nc);

tnjrxe=nan(nt,nc);
tjjrxe=nan(nt,nc);
tjuprxe=nan(nt,nc);
tjdnrxe=nan(nt,nc);

for j=1:nt
    for i=1:nc
        if indjp(j,i)<1 && resilrie(j,i)>0
            tnjlre(j,i)=resilr(j,i);
        end
        if indjp(j,i)>0 && resilrie(j,i)>0
            tjjlre(j,i)=resilr(j,i);
        end
        if indjup(j,i)>0 && resilrie(j,i)>0
            tjuplre(j,i)=abs(resilr(j,i));
        end
        if indjdn(j,i)>0 && resilrie(j,i)>0
            tjdnlre(j,i)=abs(resilr(j,i));
        end
        
        if indjp(j,i)<1 && resirxie(j,i)>0
            tnjrxe(j,i)=resirx(j,i);
        end
        if indjp(j,i)>0 && resirxie(j,i)>0
            tjjrxe(j,i)=resirx(j,i);
        end
        if indjup(j,i)>0 && resirxie(j,i)>0
            tjuprxe(j,i)=abs(resirx(j,i));
        end
        if indjdn(j,i)>0 && resirxie(j,i)>0
            tjdnrxe(j,i)=abs(resirx(j,i));
        end
    end
end
clear i j

% making hlfday window based return (CAR) and moments

hwresirxe=nan(nhlf,nc);
hwresilre=nan(nhlf,nc);
hwresivol=nan(nhlf,nc);
hwresiskew=nan(nhlf,nc);
hwresikurt=nan(nhlf,nc);
hwresijrvol=nan(nhlf,nc);
hwresijvol=nan(nhlf,nc);
hwresijupv=nan(nhlf,nc);
hwresijdnv=nan(nhlf,nc);
hwresijv=nan(nhlf,nc);

for hw=1:nhlf
    tempierx=resirxie(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempielr=resilrie(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    temprx=resirx(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    templr=resilr(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempnjlr=tnjlre(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempjjlr=tjjlre(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempjuplr=tjuplre(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempjdnlr=tjdnlre(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempnjrx=tnjrxe(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempjjrx=tjjrxe(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempjuprx=tjuprxe(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempjdnrx=tjdnrxe(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    
    for i=1:nc
        ttempierx=tempierx(:,i);
        ttempielr=tempielr(:,i);
        ttemprx=temprx(:,i);
        ttemplr=templr(:,i);
        ttempnjlr=tempnjlr(:,i);
        ttempjjlr=tempjjlr(:,i);
        ttempjuplr=tempjuplr(:,i);
        ttempjdnlr=tempjdnlr(:,i);
        ttempnjrx=tempnjrx(:,i);
        ttempjjrx=tempjjrx(:,i);
        ttempjuprx=tempjuprx(:,i);
        ttempjdnrx=tempjdnrx(:,i);
        
        ttempierx(isnan(ttempierx))=[];
        ttempielr(isnan(ttempielr))=[];
        ttemprx(isnan(ttemprx))=[];
        ttemplr(isnan(ttemplr))=[];
        ttempnjlr(isnan(ttempnjlr))=[];
        ttempjjlr(isnan(ttempjjlr))=[];
        ttempjuplr(isnan(ttempjuplr))=[];
        ttempjdnlr(isnan(ttempjdnlr))=[];
        ttempnjrx(isnan(ttempnjrx))=[];
        ttempjjrx(isnan(ttempjjrx))=[];
        ttempjuprx(isnan(ttempjuprx))=[];
        ttempjdnrx(isnan(ttempjdnrx))=[];
        
        if size(ttempierx,1)>npday/4-1 && size(ttemprx,1)>npday/4-1 && ...
                size(ttemplr,1)>npday/4-1 && size(ttempielr,1)>npday/4-1
            hwresirxe(hw,i)=mean(ttemprx)/100*npday/2;
            hwresilre(hw,i)=mean(ttemplr)*npday/2*100;
            hwresivol(hw,i)=sqrt(mean(ttemplr.^2)*npday/2);
            hwresiskew(hw,i)=(sqrt(npday/2)*mean(ttemplr.^3)*npday/2)/...
                hwresivol(hw,i)^3;
            hwresikurt(hw,i)=(npday/2*mean(ttemplr.^4)*npday/2)/hwresivol(hw,i)^4;
            
            hwresijrvol(hw,i)=sqrt(((npday/2)/size(ttempielr,1))*...
                sum(ttempnjlr.^2));
            hwresijvol(hw,i)=sqrt(sum(ttempjjlr.^2)*(npday/2)/size(ttempielr,1));
            hwresijupv(hw,i)=sqrt(sum(ttempjuplr.^2)*(npday/2)/size(ttempielr,1));
            hwresijdnv(hw,i)=sqrt(sum(ttempjdnlr.^2)*(npday/2)/size(ttempielr,1));
            hwresijv(hw,i)=hwresijupv(hw,i).^2+hwresijdnv(hw,i).^2;
            if hwresijv(hw,i)>hwresivol(hw,i)^2
                hwresijv(hw,i)=hwresivol(hw,i)^2;
            end
        end
    end
end
clear hw k i
clear tnjlre tjjlre tjuplre tjdnlre tnjrxe tjjrxe tjuprxe tjdnrxe temprxie
clear ttempielr ttempierx ttemprx ttemplr ttempnjlr ttempjjlr ttempjuplr ttempjdnlr
clear ttempnjrxe ttempjjrxe ttempjuprxe ttempjdnrxe ttempnjrx
clear tempielr tempierx temprx templr tempnjlr tempjjlr tempjuplr tempjdnlr
clear tempnjrx tempjjrx tempjuprx tempjdnrx ttempjdnrx ttempjjrx ttempjuprx

hwresivol=100*hwresivol;
hwresijrvol=100*hwresijrvol;
hwresijupv=100*hwresijupv;
hwresijdnv=100*hwresijdnv;
hwresijvol=100*sqrt(hwresijv);

% making BAS and interest rate differential
hwbas=nan(nhlf,nc);    % 12 hr BAS
tpbas=bas.*lreie;

for hw=1:nhlf
    tempie=lreie(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tempbas=tpbas(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    tpsumie=nansum(tempie,1);
    tpmeanbas=nanmean(tempbas,1);
    
    for i=1:nc
        if tpsumie(1,i)>npday/4-1
            hwbas(hw,i)=tpmeanbas(1,i);
        end
    end
    clear i tempie tempbas tpsumie tpmeanbas
end
clear hw tpbas

hwintd=nan(nhlf,nc);

for hw=1:nhlf
    tempintd=intraintd(winstart(hw,1):winstart(hw,1)+npday/2-1,:);
    hwintd(hw,:)=nanmean(tempintd,1);
end
clear hw tempintd tpbas

% using STATA change the data shape from short to long
% (i.e., on STATA, "reshape long rxe_ vol_ jrvol_ jdnvol_ jupvol_ jdnvol_ bas_ intd_ , i(num) j(fx)"
% (variables are 12 hr excess returns, total volatility, jump robust
% volatility, negative jump volatility, positive jump volatility, BAS,
% interest rate differential)
% (then read the long shape data in matlab and call "longdatasum")

%% 2. Table 1

tplong=size(longdatasum,1);

for j=1:tplong
    if temp(j,6)>1
        temp(j,6)=nan;
    end
end
clear j

sumstat(:,1)=nanmean(temp,1)';
sumstat(:,2)=nanstd(temp,1)';
sumstat(:,3)=skewness(temp,1,1)';
sumstat(:,4)=kurtosis(temp,1,1)';

corrfull=corrcoef(temp,'Rows','pairwise');

temppre=nan(tplong,7);
tempone=nan(tplong,7);
temptwo=nan(tplong,7);

for j=1:tplong
    if sum(longdatasum(j,8:10),2)>0
        temppre(j,:)=temp(j,:);
    end
    if longdatasum(j,11)>0
        tempone(j,:)=temp(j,:);
    end
    if longdatasum(j,12)>0
        temptwo(j,:)=temp(j,:);
    end
end
clear j

sumstatpre(:,1)=nanmean(temppre,1)';
sumstatpre(:,2)=nanstd(temppre,1)';
sumstatpre(:,3)=skewness(temppre,1,1)';
sumstatpre(:,4)=kurtosis(temppre,1,1)';

corrpre=corrcoef(temppre,'Rows','pairwise');

sumstatone(:,1)=nanmean(tempone,1)';
sumstatone(:,2)=nanstd(tempone,1)';
sumstatone(:,3)=skewness(tempone,1,1)';
sumstatone(:,4)=kurtosis(tempone,1,1)';

corrone=corrcoef(tempone,'Rows','pairwise');

sumstattwo(:,1)=nanmean(temptwo,1)';
sumstattwo(:,2)=nanstd(temptwo,1)';
sumstattwo(:,3)=skewness(temptwo,1,1)';
sumstattwo(:,4)=kurtosis(temptwo,1,1)';

corrtwo=corrcoef(temptwo,'Rows','pairwise');


